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ABSTRACT 

We investigate the spatial distribution of galactic satellites in high resolution simu- 
lations of structure formation in the ACDM model: the Aquarius dark matter simu- 
lations of individual halos and the Millennium II simulation of a large cosmological 
volume. To relate the simulations to observations of the Milky Way we use two alter- 
native models to populate dark halos with "visible" galaxies: a semi-analytic model 
of galaxy formation and an abundance matching technique. We find that the radial 
density profile of massive satellites roughly follows that of the dark matter halo (un- 
like the distribution of dark matter subhalos). Furthermore, our two galaxy formation 
models give results consistent with the observed profile of the 11 classical satellites 
of the Milky Way. Our simulations predict that larger, fainter samples of satellites 
should still retain this profile at least up to samples of 100 satellites. The angular dis- 
tribution of the classical satellites of the Milky Way is known to be highly anisotropic. 
Depending on the exact measure of flattening, 5-10 per cent of satellite systems in 
our simulations are as flat as the Milky Way's and this fraction does not change when 
we correct for possible obscuration of satellites by the Galactic disk. A moderate flat- 
tening of satellite systems is a general property of ACDM, best understood as the 
consequence of preferential accretion along filaments of the cosmic web. Accretion of 
a single rich group of satellites can enhance the flattening due to such anisotropic 
accretion. We verify that a typical Milky Way-mass CDM halo does not acquire its 
11 most massive satellites from a small number of rich groups. Single-group accre- 
tion becomes more likely for less massive satellites. Our model predictions should be 
testable with forthcoming studies of satellite systems in other galaxies and surveys of 
fainter satellites in the Milky Way. 

Key words: 



1 INTRODUCTION 



The satellites of the Milky Way offer a number of criti- 
cal tests of the current cosmological paradigm, the ACDM 
model. Their abundance and internal structure are sensi- 
tive to the nature of the dark matter and thus to the high 
frequency end of the linear power spectrum of density fluctu- 
ations. Their spatial distribution is sensitive to the gravita- 
tional evolution of dark matter on galactic and supergalactic 
scales. 

High resolution simulations of the formation of galac- 
tic cold dark matter halos have revealed that a large num- 
ber of substructures survive to the presen t day, account- 
ing f o r about 10% of the total halo mass C Dicmand et alj 
120071 : ISpringel et al.l 120081 . and references therein). Since 
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only about two dozen satellites are known to orbit in the 
halo of the Milky Way, this property is frequently deemed 
to pose a "satellite problem" for CDM-based cosmologies. In 
fact, it was shown about a decade ago that basic processes 
that are unavoidable during galaxy formation, such as su- 
pernova feedback and early Hydrogen reionization, readily 
explain why a visible satellite galaxy can form only in a 
tiny fraction of the surviving s ubhal os jBullock et alj|200d : 
iBenson et ail 120021 ; ISomervilld 120021 '). This result has been 
confirm ed repeatedly in recent years using semi-analytic 
models JCooper et alj|20ld: iGuo et al.ll201ll : iLi et alj|20lb1 ; 



models uuooper et ; 
IMaccio et alj l20ld : 



Font et al 



201 lib cosm o logical 



dynamical simulations (lOkamoto et al.l l20ld : iParrv et al.l 



120121 : IWadepuhl &[ Springe] 12011!) and simplified semi 



empirical models ( Kravtsov et alj|2004l : iKoposov et alj2009l : 



iBusha et alll20ld : iMufioz et alj|2009l ) 

A different kind of theoretical challenge is posed by the 
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spati al distribution of the 11 classical sat ellites of the Milky 
Way. iLvnden- Belli (119761 ). iKunkell (| 19791 ) and iLvnden-Befll 
(1 19821 ) noted that these satellites lie very close to a great 
circle on the sk y that is appro x imate ly perpendicular to the 
Galactic Plane. iKroupa et alj (|2005l ) deemed such a highly 
flattened structure to be extremely unlikely in a CDM cos- 
mology but, using N- body simula t ions of halo formation 
in the ACDM model iKanget all (120051). ILibeskind et all 
(|2005l ). IZentner et ail (|2005l ) and ILibeskind et alj (|2009l ) 
showed explicitly that this presumption is incorrect. Such 
flattened s atellite distribut i ons w ere dubbed "great pan- 
cakes" by ILibeskind et alj (|2005l ) who ascribed them to 
highly anisotropic accretion of proto-subhalos along the fil- 
aments of the cosmic web. Correlated accretion along fila- 
ments was also identified as the cause for the polar align- 
ment of satellite disks found by iDeason et all l|201lf) in 20 
percent of satellite systems (with more than 10 bright mem- 
bers each) in the "GIMIC" N-body/gasdynamic simulations 
jCrain et alj|2009l ). Although all of these studies found that 
flattening of satellite systems is common in ACDM, they 
also found that the high degree of flattening in the Milky 
Way system is atypical. 

iLi fc Helmil | 20081 ) showed that the flattening effects of 
anisotropic accretion are greatly enhanced in cases where 
infalling dark matter subhalos belong to a 'group' sharing 
similar infall times and orbital angular momentum orienta- 
tions (the members of such groups do not have to be bound 
in a single DM halo before infall). In a Mikly Way-mass sys- 
tem, they found that samples of ~ 10 subhalos drawn from 
one or two such groups readily produced configurations as 
flat as that of the classical Milky Way satellites. Extrapo- 
lating this result, they suggested that the highly correlated 
infall of a rich group of satellites was a possible explanation 
for the abnormal degree of flattening seen in the Milky Way. 
However, they could not say whether this was a probable ex- 
planation since they did not calculate the likelihood of all 
the 11 brightest satellites being members of only one or two 
such groups. 

An important limitation of the simulations that have 
been analysed so far to study the spatial distribution of 
satellites is their relatively low resolution. Low resolution 
can give rise to excessive tidal disruption and to the arti- 
ficial merging of some subhalos, potentially obscuring the 
true spatial distribution. In this paper, we use the suite of 
six simulations of individual Galactic halos from the Aquar- 
ius project, which are amongst the highest resolutio n CDM 
simulations carried out to date (|Springel et alJlioOSl ) , as well 
as a sample of similar halos from the Millennium-II simula- 
tion l|Bovlan-Kolchin et alj|2009l ). Millennium-II has lower 
resolution than Aquarius but follows halo formation in a 
cosmological volume (a cube of side 100 ft -1 Mpc, where h 
denotes the Hubble constant in units of 100 kms _1 Mpc _1 ). 
With the Aquarius simulations we are able to study satel- 
lites down to very small masses in six galactic halos and with 
the Millennium-II we are able to study the massive satellites 
of a large sample of galactic halos. Since the Aquarius ha- 
los are zoom-resimulations of regions in the Millennium-II 
volume, we are able also to study the effects of numerical 
resolution. 

We rank satellites in our simulations by stellar mass 
using two different techniques: semi-analytic modelling and 
a simple assignment of the brightest satellites to the most 



massive proto-subhalos. As we show, the two approaches 
pick out similar subsets of subhalos as satellite hosts. With 
these samples, we revisit the radial distribution of satel- 
lites and the great pancake and, for the first time, we in- 
vestigate how these properties depend on the number of 
satellites considered (in samples ranked by stellar mass). 
This allows us to make pr edictions for forth coming sur- 
veys such as Pan-STARRS (|Kaiser et alj|201oT ). which may 
discover a large population of very faint satellites in the 
Milky Way and M31, and the ongoing Pandas survey of M31 
|McConnachie et alj 120091; iMartin et alj|2009l ). We also ex- 
tend the work of ILi fc Helmil (|2008l ) by using Aquarius to 
investigate the multiplicity function of groups of massive 
satellites. 

The paper is organized as follows. In Section [21 we dis- 
cuss the suite of simulations and galaxy formation models 
that we use. In Section [3] we investigate the radial distri- 
bution of satellites, the prevalence of great pancakes, and 
groups of massive satellites. Our conclusions and discussion 
are presented in Section [4] 



2 SATELLITE IDENTIFICATION 

We begin by providing a brief introduction to the Aquarius 
and Millennium-II simulations and then describe the two al- 
ternative techniques we have employed to identify satellites 
in them. 



2.1 Aquarius and Millennium II 



The Aquarius Project (|Springel et alj 120081 ) consists of a 
suite of very high resolution N-body simulations of six dark 
matter halos of mass similar to that expected for the halo 
of the Milky Way. The simulations assume the ACDM cos- 
mology, with param eters consistent with the WMAP 1 data 
l|Spergel et alj|2003l ): matter density parameter, Qm = 0.25; 
cosmological constant term, S1a = 0.75; power spectrum 
normalization, as = 0.9; spectral slope, n 3 = 1; and Hub- 
ble parameter, h = 0.73. These values are inconsistent with 
the WMAP 7 data at about the 2a level, but this discrep- 
ancy is of no cons equence for the analysis of this paper (see 
IWang et al]|2012l ). 

The six Aquarius halos are labelled "Aq-A" through 
"Aq-F". Each was resimulated at different resolution in or- 
der to assess numerical convergence. A suffix 1 to 5 identifies 
the resolution level, with level 1 denoting the highest resolu- 
tion. In this study, we analyse primarily the level 2 simula- 
tions. For the level 2 simulations, the particle mass is m p ~ 

1 x 10 4 /i _1 Mq and the softening length is e = 48/i _1 pc. 
At z = 0, the six halos have a "virial" mass, N'hoo ~ 1- 

2 x 10 12 Mq, where M200 is the mass contained within 
?~200, the radius of a sphere of mean density 200 times the 
critical density for closure. The circular velocity curves of 
the halos peak at Knax = 220 ± 40 km s _1 . Although the 
Aquarius halos hav e similar final ma sses, they have varied 
formation histories (|Wang et alj|201ll ) . For fu rther details of 
the A quari us Project, we refer t he reader to ISpringel et alj 
(|2008h and lNavarro et alj (|2010l ). 

At every snapshot in the simulations we find nonlinear 
structures using t he friends-of-friends (FOF) algorithm of 
iDavis et alj (|l985l ). with a linking length of 0.2 times the 
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Fi gure 2. The project ed positions of satellites identified using 
the lCooper et alj BOIOT ) semi-analytic model (diamonds) and the 



^pcak method (filled circles) in the Aq-A2 halo. From top left to 
bottom right, positions are shown for the top 11, 30, 60 and 100 
satellites ranked by stellar mass in each model. The overlap ratio 
between the two samples is given in the legend. It decreases from 
/ = 0.91 for the top 11 to / = 0.7 for the top 100 satellites. 



Fig ure 1. The re l ations hip between the stellar mass assigned by 
the ICooper et al. semi-analytic model to subhalos in the 

six Aquarius simulations and the maximum value of the rota- 
tion curve attained by each subhalo over its entire past history. 
(Satellites with unresolved dark matter haloes are excluded.) The 
tight correlation for the most massive satellites motivates a simple 
^pcak model for ranking satellites by stellar mass. 



iBovlan-Kolchin et alJl2009l ). The resolution of the MRU is 
simila r to that of the simulatio ns analysed by iKang et alj 
|2005f ) and IZentner et all d2005l), but it is a f a ctor o f ~ 30 
worse than the simulation bv lLibeskind et alj (|2005h and a 



factor of ~ 20 better than that by Libeskind et al. |2009) 



mean interparticle separation and 32 particles as the mini- 
mum number of particles per group. We also identify bound 
substructu res within each FOF halo using the SUBFIND al- 
gorithm of lSpringel et al] |2005). Merger trees for su bhalos 
are constructed as described in Springel et alj (|2008l ). 

The Millennium-II run (MRU) is a cosmological sim- 
ulation in which 2160 3 particles were followed in a cu- 
bic box of side length Lbox = 100 Mpc h^ 1 . This vol- 
ume is 1 25 times smaller tha n that of the Millennium sim- 
ulation (|Springel et al]|2005l ). and the mass resolution is 
correspondingly 125 times better: each particle has mass 
6.88 x lO 6 ft _1 M0. The cosmological model is the same as 
that assumed for the Aquarius Project. For furthe r details, 
we refer the reader to lBovlan-Kolchin et alj (|2009l ). 

The Aquarius halos are resimulations of regions selected 
from a low resolution version of the MRU with identical 
large-scale density perturbations and phases in their initial 
conditions. We can therefore find the counterparts of each 
Aquarius halo in the MRU and carry out resolution tests (see 



2.2 Galaxy formation models 

There are two techniques for modelling satellite galaxies in 
detail using simulations: direct hydrodynamic simulations 
and semi-analytic modelling applied to a high-resolution N- 
body simulation of the dark matter. The two techniques 
are similar, differing primarily in the simplified, spheri- 
cally symmetric treatment of gasdynamics in the semi- 
analytic method. In this study, we are interested in very faint 
satellites. Even the highest resolutio n gasdynamic simula- 
tions of satellites performe d to date (|Okamoto et alj|20ld : 
IWadepuhl fc SpringelHioTlh can resolve only the dozen or so 
brightest galaxies pf a Milky Way system and so they are not 
suitable for our purposes. Semi-analytic models, on the other 
hand, are restricted only by the resolution of the associated 
N-body simulation so we resort instead to the semi-analytic 
model s applied to the Aqua rius simulat i ons bylCooper et alj 
(|2010h and to the MRU bv lGuo et all (120111 ). We refer the 
reader to these papers for details of the implementation. An 
added advantage of the semi-analytic models is that, unlike 
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the direct gasdynamic simulations, these models are known 
to agree with a whole range of other properties of the galaxy 
population at various epochs. 

In practice, despite the complexity of the many physi- 
cal processes they incorporate, both gasdynamic and semi- 
analytic simulations predict an approximately monotonic re- 
lation between the stellar mass of a galaxy and the total 
mass of the halo in which it forms. As we are only inter- 
ested in the stellar mass rank of satellites, it is instructive 
to consider a much simpler model for satellites in which the 
stellar mass (hence luminosity) of each satellite is assumed 
to be proportional to the highest value of the maximum of 
the circular velocity curve (Vmax) attained by the halo in 
which the satellite forms throughout its entire history. We 
denote this by Vp Ca k; it corresponds roughly to the value 
of V max just prior to the halo being accreted into the main 
halo and becoming a subhalo. We refer to this as the "Vpeak 
model". By comparing our semi-analytic results with this 
simpler but fundamentally similar model, we can check if 
our conclusions depend on the details of the semi-analytic 
treatment, which introduces a scatter in the relationship be- 
tween stellar mass and (maximum) halo mass. 

The correspondence between V D ea k and M Bta r, the ste l- 
lar mass of the semi-analytic model o f I Cooper et all (|2010h . 
is illustrated in Fig. Q] The expected monotonic relation is 
apparent. However, there is considerable scatter in the rela- 
tion for subhalos with V pC ak ^ 25kms -1 . At higher masses, 
in the range of the 11 classical satellites, the relation flattens 
and the scatter is greatly reduced. In this regime, over 90% 
of the satellites with the highest V pC ak values also have the 
highest A/star according to the semi-analytic model. Fig. Q] 
shows that our simple Vp Ca k model will give a ranking of 
subhalos by stellar mass very similar to that of the semi- 
analytic model, in most cases. 

The internal structure of halos that are close to the res- 
olution limit of an N-body simulation is, of course, poorly 
modelled. This can lead to the artificial disruption of these 
halos when they are accreted into larger halos. The baryonic 
properties of their associated galaxies can still be followed 
in semi-analytic models even after the subhalos fall below 
the resolution of the simulation, but as their orbits cannot 
be tracked within the N-body simulation their spatial distri- 
bution is uncertain and model-dependent. The resolution of 
the Aquarius simulations is so high that we need not worry 
about g alaxies with unre solved haloes (see discussion of this 
point in iFont et alJkoill ). However, at the lower resolution 
of the MRU, galaxies in unresolve d haloes do need to be 
taken into account, as discussed bv lGuo et al.l (|201ll 'l. This 
conclusion is supported by our investigation of the radial dis- 
tribution of satellites in MRU in the following section. We 
demonstrate that convergence with Aquarius is only possi- 
ble if we include galaxies without resolved subhaloes in the 
Guo et al. model. In addition, because the circular velocity 
scale corresponding to the Vp Ca k of the most massive satel- 
lites (30-60 kms -1 ) is close to the resolution limit of MRU, 
we can only make reliable comparisons between the semi- 
analytic and the Vp ea k models in the Aquarius simulations, 
and not in the MRU. 

In Fig. [2] we compare the projected positions of the top 
11, 30, 60, and 100 satellites ranked by stellar mass in Aq- 
A-2 according to the semi-analytic model (diamonds) and 
the Vp C ak method (full circles). As we can see in the top- left 
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Figure 3. The fraction of satellites ranked according to different 
criteria that overlap with the top N ga i satell ites ranked according 
to the stellar mass assigned to them by the lCooper et al] ll2010h 
semi-analytic model. The data correspond to the mean value for 
all six Aquarius halos. The black line shows the overlap with 
our standard Vpeak model and the la error bar obtained from 
the 6 Aquarius simulations. The green line shows the result of 
ranking satellites according to their present mass and the blue 
line according to the maximum mass ever attained by the halo 
in which they reside today. The red line is similar to the Vpeak 
model, except that in this case V pea k is defined as the maximum 
value attained by V max before z = 6, approximately the redshift 
of reionization. 



panel, 10 out of the top 11 satellites are the same in both 
models, an overlap ratio, / = 10/11 = 0.91. As we consider 
smaller and more numerous satellites, / decreases, but even 
for the top 100 satellites, / is still as high as 70 percent. 

In Fig. [3] we compare the overlap fractions between the 
top rigai sate llites ranked accord ing to the stellar mass as- 
signed by the lCooper et all (|2010h model and various simpli- 
fied but plausible models for identifying satellites with sub- 
halos, including our standard Vpeak model (black line with 
error bars). The overlap fraction for the Vpeak model tends 
to a constant, / ~ 0.75, for n ga i J> 60. The green line cor- 
responds to the case in which satellites are i dentified with 
the m ost massive subhalos at the present day l|Stoehr et al.l 
l2002f ). Clearly, this method fails to place satellites in similar 
subhalos to those picked out by our two standard models. 
The reason for this is simply that the present mass of a 
subhalo is significantly affected by tidal stripping and thus 
is only weakly correlated with the mass of the halo prior to 
accretion, during which time most of the satellite stars form. 

Identifying satellites instead with subhalos ranked by 
the maximum mass they ever attain (blue line) gives a better 
match to the standard models since this mass is better cor- 
related with the satellite's stellar mass or luminosity. Even 
in this case, though, the overlap falls to just over 60 percent 



© 2009 RAS, MNRAS 000. \TWL2\ 



Spatial Distribution of Galactic Satellites 5 




50 100 150 200 250 
r [kpc] 




50 



1 00 150 
r [kpc] 



200 250 



Figure 4. The radial profiles of the top If satellites in differ- 
ent models. Results have been averaged over the 6 Aquarius ha- 
los. The profile of the ff classical satellites of the Milky Way is 
shown by the filled circles joined with a thick dashed black curve. 
The radial profile of the dark matter is shown by the black dot- 
ted curve and the radial profile of all resolved subhalos in the 
six Aquarius halos by the red dotted curve. The green, blue, red 
and black lines show the distribution of satellites in the Aquarius 
simulations, identified using different models, as indicated in the 
legend. The green diamonds with error bars correspond to satel- 
lites in the semi-analytic model applied to the MRU, including 
satellites without resolved dark matter haloes. 



for large satellite populations. The red line shows what hap- 
pens if we identify satellites according to the value of Vpeak 
attained prior to z = 6, roughly the redshift of reionization. 
This model gives similar results to our standard Vpeak model, 
except for the most massive satellites, reflecting, in part, the 
fact that reionization plays a relatively mino r role in setting 
the stellar mass of th ese galaxies (see e.g. lOkamoto et al.l 
l20ld ; lFont et al.|[201lh . We conclude that our choice of Vpeak 
is the most appropriate of these alternatives and in mos of 
the remainder of this paper we focus on comparisons be- 
tween satellites ranked by this property and by the stellar 
mass in our semi-analytic models. 



3 RESULTS 

In this section we consider the spatial distribution of satel- 
lite galaxies predicted in both our semi-analytic and Vp ea k 
models and compare them with data for satellites in the 
Milky Way. 



3.1 Radial profiles 

The radial distribution of the II classical satellites of the 
Milky Way is shown by the filled circles in Fig. [4] The Galac- 



Figure 5. The radial profile of the II brightest satellites accord- 
ing to the Vpeak model in each of the six Aquarius halos (solid 
lines) and in their counterparts in the MRU (dotted lines). The di- 
amonds and large error bars correspond to the average profile and 
lcr dispersion of the top If satellites identified by the Vpeak model 
in f686 MRU halos of mass 1 X lO 12 Af < M 20 o < 2 X IO 12 M . 
The profile of the If classical Milky Way satellites is shown by 
the black circles joined with a dashed black curve. 



tocentric distance used here is taken from iMated (|l99Sl ). 
It is close to the radial distribution of the dark matter in 
the Aquarius halos, shown in the figur e by the black dot - 
ted line, a similarity a lread y noted by Kang et al.l (|2005T ) . 
iLibeskind et all |2005f ) and Font et all (|201lh . To compare 
the data to our models, we select the top 11 satellites in each 
model that lie within a distance of 250 kpc from the centre 
of the halo, roughly the distance of Leo I, the most distant 
of the observed Milky Way eleven. 

Results for the semi-analytic, Vpeak, Vpeak {z > 6) and 
M m ax(z = 0) models are shown by the green, black, blue 
and red solid lines respectively. All, except the last, are 
roughly consistent with the data within the la rge uncertain- 
ties due to small num ber statistics (see also ICooper et al.l 
|2010| ; iFont et all l201ll ). By contrast, identifying satellites 
with the most massive present day subhalos (red line) leads 
to a much more extended distribution which, as first noted 
bv lLibeskind et all (|2005l ) , closely follows the radial profile of 
th e entire halo populatio n. This is consistent with the finding 
by ISpringel et al.l (|2008l ) that the radial profile of subhalos 
depends very weakly on subhalo mass. The other, more real- 
istic, models place bright satellites in subhalos whose radial 
profiles are strongly biased relative to that of the subhalo 
population as a whole and which happen to lie close to the 
dark matter radial profile. We also show results for Milky 
Way galaxies from the semi-analytic model applied to the 
MRU, taking care to include satellites which survive even 
though their dark matter halos have been tidally stripped. 



© 2009 RAS, MNRAS 000.ITHT21 



6 



Jie Wang, Carlos S. Frenk and Andrew P. Cooper 



v 



c 

o 



u 
o 



c 
o 

A- 1 

u 
o 




50 100 150 200 250 50 100 150 200 250 
r [kpc] r [kpc] 

Figure 6. The cumulative radial profiles of the top -/V ga i satellites 
in each of the Aquarius halos in the semi-analytic (solid lines) 
and Vpeak (dotted lines) models. The profile of the 11 classical 
satellites of the Milky Way is shown in all panels by the solid 
dots joined by a thick dashed black curve. 



The distrobution is close to the observational data. However, 
had we not included satellites without resolved dark matter 
halos, the profile would lie close to that of all subhaloes (red 
dotted). 

The model lines in Fig. [4] are smoother than the Milky 
Way data because we averaged the distributions of the 11 
top-ranked satellites from each of the 6 Aquarius halos. 
There is, in fact, considerable halo to halo scatter, as may 
be seen in Fig. O which shows the radial profiles of the 11 
brightest satellites predicted by the semi-analytic model in 
each halo. Five of these lie very close to the Milky Way data 
but in Aq-A, the radius enclosing half the satellites is over 
twice as large as in the other 5 Aquarius halos. 

Fig. [S] also highlights the importance of high resolu- 
tion simulations for this kind of study. Since the Aquar- 
ius halos were selected from the MRU, we can iden- 
tify their counterpart s in th at simulation, as described by 
iBovlan-Kolchin et alj (|2010l ). The radial profiles of the top 
11 satellites in these counterparts, as identified in the Vp ea k 
model, are shown by the dotted lines in the figure, with 
colours corresponding to those of the Aquarius halos. Satel- 
lites in the MRU counterparts are clearly much more ex- 
tended than the Aquarius subhalos. This is typical of the 
satellite distributions in Milky Way-mass halos in the MRU, 
as illustrated by the diamonds with error bars, obtained by 
averaging the distributions of the 11 brightest satellites in 
the 1686 MRU halos of mass 1 x 10 12 < M 2 oo < 2 x 10 12 . 
The MRU has about 1000 times lower resolution than the 



level 2 Aquarius simulations. As a result, subhalo orbits are 
not followed as accurately and their disruption timescale is 
artificially reduced. 

The semi-analytic model of iGuo et all (|201ll ) allows 
galaxies to survive after their host subhalo falls below the 
resolution limit. A position is assigned to these galaxies by 
tracking the most-bound particle of the host subhalo from 
the time it was last resolved. This position is unlikely to be 
a very accurate estimate of the true orbit of the satellite; 
for example, Guo et al. make an analytic correction to this 
single-particle orbit when calculating the tidal disruption of 
the satellite, as the orbits of individual particles do not decay 
through dynamical friction in an N-body simulation. Never- 
theless, when we include satellites with unresolved haloes in 
our top 11 sample, as was done in Fig. 2] (but not in Fig. O, 
and use their single-particle position estimates, the MRU 
radial profiles for the Aquarius galaxies are in reasonable 
agreement with their higher resolution counterparts. 

As we expect that larger and more complete samples 
of fainter satellites will soon be available for the Milky Way 
and other galaxies, it is interesting to know how the radial 
profiles of satellites are predicted to vary with the luminos- 
ity cut of the sample. This is a question that we can address 
with the Aquarius simulations and we do this in Fig. [(J 
where we show the radial profiles of the top 11, 30, 60 and 
100 galaxies predicted in the semi-analytic (solid curves) 
and Vp C ak (dotted curves) models. The profiles depend only 
weakly on the size of the sample within the mass range we 
can examine: they are all similar to that of the 11 classical 
satellites of the Milky Way. The exception is Aq-A whose 
profiles become increasingly less discrepant with the other 
5 Aquarius cases as the satellite sample size increases (al- 
though it is still the least concentrated even when 100 satel- 
lites are considered). By most other m easures Aq-A is repre- 
senta tive of haloes of the same mass (|Bovlan-Kolchin et al.l 
2010). suggesting that the top 11 satellites may often exhibit 
'unusual' configurations in otherwise 'typical' haloes. 



3.2 The angular distribution of satellites and the 
great pancake 

As we have seen in the preceding section, the semi- 
analytic and Vpeak models reproduce the observed ra- 
dial profile of the 11 classical satellites of the Milky 
Way. In this subsection we investigate if they also repro- 
duce the highly flatt ened configu r ation seen in the Milky 
Way, first noted by Lvnden-Belll ( 1976 ) and investigated 



more r ecently by 



J2005I); iKane et al 
iDeason et alj feoil 



Kroupa et al.l ( 2001 ); iLibeskind et al] 



( 20051 ); ILibeskind et all (|2009l ) and 



We use simila r tec hniques to those in - 



troduced by ILibeskind et all (|2005l ) and iKang et all (|2005l) 

to describe the flattening of the satellite population, namely, 
the overall shape of its distribution and its thickness in the 
flattest dimension. 



3.2.1 The axis ratios of satellite systems 

To determine the shape of a satellite system we calculate the 
moment of inertia tensor of its members (weighting each of 
them equally) and derive the principal axes of the distribu- 
tion. The eigenvalues of the diagonalised inertia tensor are 
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Figure 7. The axis ratios describing the shape of satellite systems 
in the Aquarius simulations and in the Milky Way. c/a is the 
minor to major axis ratio and b/a the intermediate to major axis 
ratio. Results for the top 11 satellites in the semi-analytic and 
^peak models are shown by asterisks and filled circles respectively, 
with error bars spanning the range obtained for a large sample of 
satellites taking into account systematic effects introduced by the 
zone of avoidance (see text for details). The axis ratios for the 11 
classical satellites of the Milky Way are shown with a red cross 
with the value of c/a highlighted by the horizontal dotted line. 
The average axis ratios of the 6 Aquarius dark matter halo (within 
f"20o) are shown as the large open triangle towards the top-right 
corner. The red square shows the effect of including Canes Ventici 
in the observational sample. 



proportional to the rms deviation of the x, y and z coordi- 
nates of the system. Denoting the major, intermediate and 
minor axes by a, b and c respectively (a > b > c), the flat- 
tening of the system may be quantified by the axial ratios 
c/a and b/a. 

The axial ratios of the system consisting of the 11 clas- 
sical satellites of the Milky Way are indicated by a red cross 
in Figs. [7] and [8] with the value of c/a — 0.18 highlighted 
by the horizontal dotted line. This i s sligh tly lower than the 
value obtained bv lLibeskind et all (|2005l ) who had an error 
in their computer code (N. Libeskind, private communica- 
tion). The faintest of the 11 classical satellites are Draco and 
Ursa Minor (My — —8.8), but Canes Venitici is just 0.2 mag 
fainter than this (|McConnachiell2012r i. We have checked that 
including Canes Venitici has virtually no effect in the esti- 
mated axis ratios of the system, as shown by the red square 
in Fig. [3 

The corrected value of c/a for the Milky Way's classical 
satellite appears quite extreme. However, with a sample of 
only 11 objects, it is not immediately obvious that a value 
even as low as c/a = 0.18 is sig nificant. To asses s the s ignif- 
icance of this result, we follow iLibeskind et~al] (|2005l ) and 
create a large number of artificial systems of 11 satellites 



Figure 8. The distribution of shape parameters for systems of 11 
satellites. The black lines show results for systems around 1686 
Milky Way like halos in the MRU and the dotted lines indicate the 
results after correcting for missing satellites in the galactic plane 
as described in the text. The red lines show results for 10000 
randomly selected systems with the same radial distribution as 
the observed satellites of the Milky Way, but with a uniform dis- 
tribution in solid angle. The three contour levels correspond to 
fractions of 30, 60 and 90 percent. The axis ratios of the 11 clas- 
sical satellites of the Milky Way are indicated by the red cross, 
with the value of c/a highlighted by the horizontal dotted line. 



in which the radial distances of the members are the same 
as the distances of the real satellites, but the latitude and 
longitude are assigned at random on the surface of a sphere. 
Isocontours of axial ratios for a set of 10,000 such isotropic 
systems, corresponding to number fractions of 30, 60 and 
90 percent, are shown as red curves in Fig. [8] We find that 
only 98 of these samples result in values of c/a smaller than 
0.18. Thus, we conclude that the distribution of the 11 clas- 
sical satellites of the Milky Way is highly anisotropic, with 
only a 1 per cent probability that such a low value of c/a 
would result from a statistical fluctuation of an intrinsically 
isotropic system. 

Applying the same procedure to the simulations we ob- 
tain the axis ratios for systems consisting of the top 11 satel- 
lites, ranked by stellar mass, that lie within 250 kpc of the 
centre in our two models. The results are shown in Fig. [7] for 
the six Aquarius halos, using asterisks for the semi-analytic 
model and circles for the Vpeak model. We also plot the av- 
erage axis ratios of the 6 Aquarius halos as a whole (within 
T20o) as a large triangle. The axis ratios of the 11 classical 
satellites of the Milky Way are plotted as a red cross. Three 
conclu sions emerge from this figu re. Firstly, as note d ear- 
lier bv lLibeskind et all (|2005l) and lKang et al.l (|2005l ). mas- 
sive satellite systems can be much flatter than the halo as 
a whole. Secondly, the two methods for ranking galaxies by 
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Figure 9. Shape parameters for systems consisting of the top 
iVg a i satellites according to our two models of satellite galaxy 
formation. The symbols (asterisks for the semi-analytic model 
and filled circles for the Vpeak model) give the directly measured 
mean values, averaged over the 6 Aquarius halos. The error bars 
span the minimum and maximum values amongst the 6 halos 
after correction for satellites missing due to obscuration in the 
Galactic plane. The four panels show results for iV ga i = 11, 22, 44, 
and 66. In the top-left panel, the values for the top 11 satellites 
of the Milky Way ranked by luminosity are indicated by a red 
cross, with the value of c/a highlighted by a red horizontal solid 
line. The triangle near the top-right corner of each panel shows 
the mean axis ratios of the entire dark matter distribution in the 
6 Aquarius halos. 



stellar mass give very similar results. Thirdly, the satellite 
systems in the simulations appear generally less flattened 
than the 11 classical Milky Way satellites. 

Although it is often assumed that the 11 classical satel- 
lites represent a complete sample of bright satellites in the 
central regions of the Milky Way, it is, of course, possible 
that bright satellites might remain undetected in the zone 
of avoidance as a result both of large extinction and confu- 
sion by foreground stars. Assuming that satellites are intrin- 
sically uniformly distributed in (1 — sin |6|), where b denotes 
Galactic latitude, IWillman et alj (|2004f ) estimated that the 
known population of dwarf satellites could be 33 percent 
incom plete because of obscuration. However, IWhiting et al.l 
argued that, at most, there are one or two dwarfs still 
undiscovered near the Galactic plane. 

To assess the possible effects of obscuration, we assume 
a geometry for the Galactic disk in our simulations and ex- 
clude satellites in the corresponding zone of avoidance. The 
"great pancake" in the Milky Way lies perpendicular to the 
Galactic disk. We therefore imagine that each of our Aquar- 
ius halos has a "Galactic plane" whose normal is parallel to 
the major axis of its satellite system. We then calculate the 
latitude of each satellite relative to this plane. If the latitude 



is less than a critical value, ±# cr it, we regard this satellite 
as undetected and exclude it from our sample of 11, replac- 
ing it by the next in the list. (For a 33 percent occulted 
sky fraction, # cr it = 9.5°.) We repeat this process iteratively 
2000 times, each time recalculating the shape of the current 
satellite system and the position of the obscuring "Galac- 
tic plane." In 2000 cases, satellites were hidden behind the 
model galactic plane and replaced, on average, with 2.6 new 
satellites. 

The results of this procedure, for 6> cr it = 9.5°, are re- 
flected in Fig. [7] as error bars spanning the range of c/a, 
with the uncorrected sample shown as a filled symbol. It is 
striking that the inferred shape of the satellite system is ex- 
tremely sensitive to the inclusion or exclusion of only one or 
two members. Allowing for sample incompleteness improves 
the comparison between models and data noticeably, with 
several of the models now coming close to the data. 

To improve our statistics, we also analysed satellite sys- 
tems in the MRU. As discussed earlier, the V poa k model 
is strongly affected by the relatively low resolution of the 
MRU, but not so the semi-analytic model, provided satellites 
with dark matter haloes below the resolution limit are in- 
cluded in the analysis. As before, we consider all dark matter 
halos with virial mass 1 x 10 12 ^ M 2 qg «S 2 x 10 12 that have 
at least 11 satellites within 250 kpc. This provides a sample 
of 1686 Milky Way halos. The axis ratios of the brightest 11 
satellites in each of the systems are shown as the black con- 
tours in Fig. [8] where levels correspond to number fractions 
of 30, 60 and 90 percent. We find only 101 satellite systems 
with an axis ratio c/a smaller than the observed value of 
0.18, giving a probability of finding a configuration flatter 
than that in the MW of only 101/1686 = 6%. 

Just as for the Aquarius sample, we make a correction 
for the possible omission of satellites near the Galactic disk. 
For each halo, we assume a Galactic plane perpendicular 
to the minor axis of the satellite distribution, draw 2000 
random samples (replacing lost satellites and iterating) and 
obtain the mean values of the axis ratios. As the dotted 
contour in Fig [8] shows, the probability of finding a system 
in the simulations as flat or flatter than that of the Milky 
Way is not sensitive to our corrections for obscuration by the 
Galactic disk (under the assumption that the minor axes of 
the satellite system and the disk are perpendicular). 

The high mass resolution of the Aquarius simulations 
allows us to investigate how the flattening of the satellite 
system varies as increasingly faint satellites are included 
in the sample and thus to make predictions for the larger 
samples that may become available in the future as sur- 
veys such as Pan-STARRSl and LSST discover new, fainter 
satellites in the Milky Way. The axial ratios for systems of 
-Wgai = 11, 22, 44 and 66 satellites are shown in the four pan- 
els of Fig. [9] The minimum and maximum values amongst 
the 6 Aquarius halos, after correcting for the effects of the 
zone of avoidance, as discussed above, are indicated by the 
error bars. Fig. [9] shows that, as the size of the satellite sam- 
ple increases, the distribution becomes i ncreasingly l e ss flat - 
tened, in agreement with the results of [Kang et al. | (|2005h . 
When the sample size reaches 66, the shape of the satellite 
configuration is very close to that of halo dark matter, indi- 
cated by a large triangle in each panel. This prediction from 
our simulations should be readily testable with potentially 
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Figure 10. The flattening (defined as the ratio of the rms height 
of the satellites relative to the best fit plane divided by the radial 
extent of the system) of the satellite distribution in the six Aquar- 
ius halos. Different panels correspond to samples containing the 
top iVg a i satellites according to our two models of satellite forma- 
tion, which are illustrated by the filled symbols, as indicated in 
the legend. The open symbols show the average of 20,000 artificial 
samples each with the same radial distribution as the correspond- 
ing Aquarius sample but with random angular positions on the 
sky. In the iV ga i = 11 panel the flattening of the 11 brightest 
satellites of the Milky Way, 0.074, is indicated by a red solid line. 
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Figure 11. The distribution of Battenings (obtained from the 
rms height of satellites relative to the best-fit plane) for satellite 
systems. The solid black line shows results for the 1686 halos in 
the MRU with 1 X 10 12 < M 2 00 =S 2 X 10 12 that have at least 
11 satellites within 250 kpc. The red histogram corresponds to 
10000 samples with the same radial distribution as the 11 classical 
satellites of the Milky Way and a random distribution of angles 
on the sphere. The dotted histogram shows results for the MRU 
corrected for the possible loss of satellites near the Galactic disk 
taken from 2000 realizations. The flattening of the 11 classical 
satellites of the Milky Way, 0.074, is shown as a dashed vertical 
line. 



forthcoming samples of new faint satellites, provided their 
selection functions are well understood. 



3.2.2 The thickness of the satellite "disk" 

An alternative method to estimate the flattening of the 
"Great pancake" is to measure t he thickness of t he sla b de- 
fined by th e panc ake. Following iKroupa et al.l (|2005l ) and 
iKang et al.l (|2005l ). we find the best fit plane to a given 
satellite sample by minimizing the root-mean-square of the 
height of each satellite relative to the plane. The thickness 
of the slab is taken to be the rms height about the best- 
fit plane and the ratio of the thickness to the radial extent 
of the system, R C ut = 250 kpc, is used to characterise its 
flattening. 

The flattening of the satellite disk in each of the six 
Aquarius halos is shown in Fig. 1101 The four panels corre- 
spond to systems with different numbers of members, N gg _i. 
Results for our two different satellite galaxy formation mod- 
els are indicated by the filled symbols. The open symbols 
show results for 20,000 artificial samples with the same ra- 
dial distribution as the corresponding Aquarius sample but 
with random angular positions on the sky. In some cases (e.g. 
Aq-A in all panels, Aq-D in the JV ga i = 11 panel, Aq-C in the 
ATgai = 22 panel), the simulated satellite systems are signif- 



icantly anisotropic but in others, the open symbols lie close 
to the filled symbols, indicating that these systems show no 
significant flattening according to this test. In the top left 
panel, we compare the model results with the measurement 
for the 11 classical satellites of the Milky Way, whose thick- 
ness (indicated by the red line) is 0.074. This population 
is clearly flatter than most of the simulations, although the 
satellite systems in Aq-C and Aq-D come close. The remain- 
ing panels in Fig.[l0]show that, as was the case for the axial 
ratios, systems containing more and more satellites become 
progressively less flattened. 

We also estimated the thickness of satellite systems in 
the 1686 Milky Way like halos found in the MRU as de- 
scribed in Section f3. 2. II which provides better statistics than 
the Aquarius halos. The histogram of flattening values is 
shown in Fig. 1111 The value measured for the 11 classical 
satellites of the Milky Way, 0.074, is indicated by the vertical 
dashed line. This degree of flattening is highly significant: of 
10000 artificial isotropic samples of 11 satellites having the 
same radial distribution as the observed satellites (red his- 
togram), only 107 are flatter than the Milky Way system. 
This corresponds to a probability of only 1.07 percent. In 
the MRU a total of 219 cases have flatter satellite systems 
than the Milky Way, corresponding to a probability of 13 
percent. When the effects of the Galactic disk are taken into 
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Figure 12. Histogram of multiple accretion satellite events in 
the six Aquarius halos. Satellites identified by the semi-analytic 
model are shown as diamonds, while those identified by the Vp ea k 
model arc shown as crosses. In both cases, red symbols corre- 
spond to the top 11 satellites and black symbols to the top 60. 
The fractions of correlated halos relative to the total number of 
satellites, / c ,sa(^V), for the semi-analytic model, and / c ,vp(-?V) 
for the Vp ea k model, where N is the number of satellites in the 
system, are given in the legend. 



account, the distribution of flattenings is not significantly 
affected (dotted histogram in Fig. 

3.3 The multiplicity of accreted groups 

As we have seen, the distributio n of satellites in the M ilky 
Way is significantly anisotropic. iLibeskind et alj |2005l ) as- 
cribed this striking property to the accretion of satellites 
along a few filaments of the cosmic web. However, in agree- 
ment with previous studies, we have shown that, although 
flattened satellite systems are common, the Milky Way's 
satellites are even m ore flattened than the expectation for 
halos of similar mass. iLi fc Helmil (|2008l ) noted that the nat- 
ural flattening induced by this anisotropic infall is strongest 
when considering samples of galaxies accreted in 'groups', 
defined by common infall conditions (time, direction and 
orientation of orbital angular momentum). Considering all 
dark matter subhalos regardless of their mass, they found 
that samples of 11 subhalos drawn from only one such in- 
falling 'group' have a configuration as flat as the classical 
Milky Way satellite distribution in 73 percent of cases. This 
probability falls to only 20 percent if only 5/11 satellites 
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Figure 13. The fraction of the top -/V ga i satellites which are ac- 
creted into the main halo in a group with at least one partner, 
as a function of -/V ga i. The dotted lines show results for satellites 
in the semi-analytic model (/ c sa)i th e solid lines for satellites in 
the Vp eak model (/ Ci y P ). 



are drawn from a single group. However, since Li & Helmi 
studied only one high-resolution halo, they were unable to 
quantify how likely it is that the most massive 11 satellites 
of the Milky Way fell in in just one or two groups and so to 
test whether or not this process is a probable explanation 

for their unusually flat configuration. 

In order to explore the effect studied by ILi fc Helmil 

(|2008f ) in our six 'typical' Milky Way halos, we identify satel- 
lites at the final time in the simulation and trace them back 
to the friends-of-friends group to which they belonged be- 
fore they were accreted. If two or more satellites were part of 
the same group at one or more subseq uent simulation sn ap- 
shots, we regard them as correlated (|Li fc Helmil |2008| did 
not require their correlated subhalos to be members of the 
same FOF group, only to share similar infall times, direc- 
tions and orbital angular momentum orientations; in prac- 
tice the two definitions are likely to have similar results). 
In Fig. [12] we plot, for each Aquarius halo, the number of 
groups in which a given number of satellites were brought 
into the main halo (counting only satellites surviving to red- 
shift zero). As before, we identify satellites using both our 
semi-analytic and Vpcak models (indicated by diamonds and 
crosses respectively) and we consider samples of A^ai = 11 
(red symbols) and A^gai = 60 (black symbols) satellites. The 
fractions of correlated halos relative to the total number of 
satellites (/ c ,sa for the semi-analytic model and / c ,vp for 
the Vpcak model) are given in the legend. 

Focussing first on the top 11 satellites in each model 
(red symbols), we can see that, except in the semi-analytic 
model for halo F, fewer than half of these satellites are mem- 
bers of groups according to our definition; that is, more than 
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half are accreted singly into the main halo. This suggests 
that the probability of the top 11 satellites falling into the 
halo as members of only two or three groups is low; the mod- 
erate Battenings of the Aquarius satellite systems represent 
the typical outcome of anisotropic accretion, and are not 
enhanced by closer orbital associations between particular 
satellites. Testing whether or not the extreme flattening of 
the Milky Way's system can be explained by the rare occur- 
rence of a ha ndful of 'rich ' or 'to p heavy' accretion events, as 
proposed bvlLi fc Helmil | 20081 ) . will require a much larger 
suite of simulations comparable to Aquarius. 

The prevalence of strongly correlated accretion in- 
creases as we consider larger satellite samples. For example, 
in halos B and F more than 15 of the top 60 satellites were 
part of a single group before accretion, as may be seen in 
Fig. 1121 The increase in / c ,sa and / c ,vp with the size of the 
satellite sample is illustrated in Fig. Qj|] and is seen in both 
our methods for ranking satellites. For systems of only 10 
satellites, the correlated fraction is ~ 10 — 50 per cent, but 
this rises to ~ 30 — 80 per cent for syste ms of 100 satellites . 
Our findings reinforce the conclusion of lLovell et al.l (|201ll ) 
that a more complete census of satellite motions in the Milky 
Way and M31 will reveal a significant fraction with common 
orbital planes. 



4 CONCLUSION AND DISCUSSION 

The spatial distribution of satellites in the Milky Way and 
other galaxies reflects both the nature of the dark matter 
and the processes of galaxy formation. In this paper we have 
used high mass resolution cosmological simulations of struc- 
ture formation in the ACDM model: the Aquarius resimu- 
lations of galactic scale dark matter halos and the Millen- 
nium II simulation (MRU) of a cubic volume 100 /i -1 Mpc 
on a side, with about 2000 halos of mass comparable to the 
Milky Way's halo. We have further employed two models of 
galaxy formation to tr ace satellites in the simu lation: the 
semi- a nalytic models of lCooper et al.l (|2010h and lGuo et al.l 
|201ll ) for Aquarius and MRU respectively and a less sophis- 
ticated model in which galaxy stellar mass is assumed to be 
proportional to Vpeak, the highest circular velocity attained 
by the halo throughout its life. Our two methods place 10 of 
the most massive 11 satellites in the same Aquarius halos. 

The combination of high-resolution simulations and ro- 
bust galaxy formation models allows us to identify satel- 
lites in a reliable way. We find that good mass resolu- 
tion is essential to obtain an accurate estimate of the ra- 
dial distribution of satellites. Indeed, even at the resolu- 
tion of the MRU, m p = 6.9 x 1O 6 /i _1 M0, many genuine 
substructures are artificially destroyed, as may be seen in 
Fig. [5] In this case, the Vpeak method for assigning stellar 
mass to satellites, as well as all previous c osmological sim- 
ulations which have po orer resolution (e.g. iLibeskind et al.l 
120051 : iKang et aTll2005l ). would give inaccurate results. By 
contrast, our semi-analytic model, which continues to track 
satellites even after they have lost their halos does give 
a faithful prediction of the theoretical expectations in the 
ACDM model. In the case of the Aquarius simulations, the 
resolution is good enough tha t galaxies with un resolved ha- 
los galaxies are unimportant (|Font et al.ll201ll ) and we can 



compare our semi-analytic results with the V pea .k method. 
The results from the two ranking methods agree well. 

Our main results concerning the satellite radial distri- 
bution are: (i) there is substantial halo-to-halo scatter in 
the radial distribution of satellites in the Aquarius simula- 
tions, but in all cases the massive satellites have a similar 
radial density profile to that of the dark matter; (ii) both the 
semi-analytic and Vp Ca k models give results in good agree- 
ment with the measured radial distribution of the 11 classi- 
cal satellites of the Milky Way; (iii) the radial density pro- 
file of larger satellite samples, going down to lower masses, 
is similar to that of the 11 most massive, at least up to 
samples of 100 satellites. This prediction may be tested by 
future surveys such as Pan-STARRS. 

We also investigated the angular distribution of satel- 
lites. This is interesting because the Milky Way's 11 bright- 
est satellites show a very aniso t ropic distribution, the "Great 
pancake" of ILibeskind et alj (|2005h . We characterized the 
angular distribution by measuring both the axial ratios of 
the satellite system and the thickness of the slab in which 
the satellites are concentrated. We verified that the pecu- 
liar distribution seen in the Milky Way is highly significant: 
only 1 percent of isotropic samples with the same radial dis- 
tribution are flatter than the Milky Way's system. In the 
1686 halos in the MRU which have a mass plausibly simi- 
lar to that of the Milky Way we find systems flatter than 
that of the classical satellites with probabilities of 6 and 
13 percent according to the axial ratio and slab thickness 
tests respectiv ely. These probab i lities are slightly lower than 
that found by ILibeskind et ahl (|2005T l. who used a slightly 
incorrect value for the flattening of the Milky Way satel- 
lites. Our probabi lities are also lower than those found by by 
IKang et alj (|2005r ) , a discrepancy that is readily understood 
as a consequence of the poor mass resolution of the simula- 
tions used by these authors. The Milky Way system thus ap- 
pears to be in the tail of the flattening distribution expected 
for massive satellites in the AC DM model. The presence of a 
satellite as massive as the LMC jBovlan-Kolchin et alj|201(il: 



iBusha et al.ll201ll : iGuo et al.ll2oTll:lLiu et al 



polar alignment of the system (|Deason et al 



2011 



and the 
2011T ) seem to 



be comparably rare, but still consistent with the predicted 
distributions. 

Finally, we investigated the extent to which satellites 
are accreted in groups. This is interesting for several rea- 
sons, including the possibility that the Great pa ncake might 
be ex plained by multiple accretion in 2-3 groups (|Li fc Helmil 
120081 ). Our simulations confirm that this is a rare occurrence 
in ACDM. On average, only 30 percent of the top 11 satel- 
lites in the Aquarius simulations share a friends-of-friends 
group with another top 11 satellite before infall; the rest 
come into the main halo without any companions of com- 
parable mass. However, multiple accretion becomes increas- 
ingly important for larger, fainter samples of satellites. For 
example, in samples of the 60 most massive satellites in two 
of the Aquarius halos (Aq-B and Aq-F), we find that as 
many as 20 come into the main halo in a single group and 
as many as 11 in the other simulated halos. This interesting 
prop erty may potent ially have observational consequences 



(e.g. ISales et aflboill ') 



ILibeskind et alj (120051 s ) proposed that filamentary ac- 
cretion is responsible fo r the flatten i ng o f Milky Way- 
like satellite populations. lLovell et al l (|201ll ) subsequently 
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showed that the subhalos of the Aquarius simulations have 
strongly correlated orbital angular mome nta as the result 
of anis otropic accretion. Recently, however. IPawlowski et all 
(|2012h re-analysed the angular moment um orientations o f 
subhalos in Aquarius (as presented in lLovell et al.l 1201 if ) 
and concluded that anisotropic accretion is unimportant for 
producing flattene d satellite system s ; in di rect c ontradiction 
with t he results of iLibeskind et aL I l|2005l ) and lLovell et al.l 
i|201lT ). They required that the Milky Way flattening should 
be reproduced in the mean, rather than being merely consis- 
tent with the expected distribution. As the Milky Way is in 
the ~ 10 per cent tail of the distribution predicted by ACDM 
(as our study and others have shown), they reject a CDM 
model in favour of tidal galaxy formation, which they claim 
readily produces highly correlated di sks of satellites. We be - 
lieve that the explanation offered by ILibeskind et all (|2005r i 
remains the most appropriate. A single randomly chosen sys- 
tem such as the Milky Way is extremely unlikely perfectly 
to represent the mean value of every measurable property. 
A larger sample of satellites aro und other galaxies will t est 
the tidal formation hypothesis of IPawlowski et ail (|2012T l in 
which highly flattened configurations are easily achieved and 
should therefore be the norm. If, on the other hand, the 
CDM model is a realistic description of nature, then the 
average satellite configurations should be only moderately 
flattened, as illustrated in Figs. [7] and [8] 

Our simulations are useful not only to test our models 
against data for the classical satellites of the Milky Way, as 
we have done here, but also to make predictions for future 
surveys that will quantify the spatial distribution of larger 
and fainter satellite samples, both in the Milky Way and in 
other galaxies. For a sample of Milky Way satellites complete 
to very fa int magnitudes, w e have shown that the luminosity 
function (|Font et al.l 1201 ll ^ . radial distribution and system 
shape should vary less from halo to halo than is the case 
for the most massive eleven. Measuring these three highly 
characteristic properties of the satellite system test both the 
ACDM cosmology and models of galaxy formation in novel, 
interesting ways. 
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